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Abstract 

' We discuss stationary solutions of the nonlinear Schrodinger equation 

(NSE) applicable to several quantum spin, electron and classical lattice sys- 
^ ■ terns. We show that there may arise chaotic spatial structures in the form of 

incommensurate or irregular quantum states and trajectories in space. 

As a first (typical) example we consider a single electron which is strongly 
coupled with phonons on a ID chain of atoms. In the adiabatic approxima- 
tion the system is conventionally described by a discrete set of NSEs. Another 
apt example is that of superconducting states in layered superconductors de- 
scribed by the same NSE. Amongst many other applications the typical ex- 
ample for a classical lattice is a system of coupled nonlinear oscillators. 

We reformulate this discrete NSE to the form of a 2D mapping. By this 
we may investigate a quantum problem by methods conventionally applied 
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to classical chaotic dynamics. We find three types of solutions: periodic, 
quasiperiodic and chaotic. We then develop a procedure which allows us to 
obtain numerical solutions of the NSE directly. This procedure may be used 
to any arbitrary accuracy and so these solutions are exact to the degree of 
precision specified. Both methods give a consistent result. When applied to 
our typical example we find that the wave function of an electron on a de- 
formable lattice (and other quantum or classical discrete systems) may exhibit 
incommensurate and irregular structures. 

PACS numbers: 02.60. Cb, 63.20.Kr, 73.20.Dx, 74.80.Dm 

1 Introduction 

Chaos is an important branch of nonlinear dynamics because chaotic behaviour 
seems to be universal It is present in mechanical oscillators, electrical circuits, 
lasers, nonlinear optical systems, chemical reactions, nerve cells, heated fluids and 
weather systems. Even more importantly, this chaotic behaviour shows qualita- 
tive and quantitative universal features which are independent of the details of the 
particular system. 

In Quantum Mechanics instead of initial conditions we have boundary conditions. 
It is commonly believed that in quantum systems chaotic structures cannot arise. 
However, the question does arise if it would be possible to have in quantum systems 
a situation analogous to classical chaotic structures? In other words, if it would be 
possible for two slightly different boundary conditions or some physical parameters 



2 



(for example, coupling constant) to correspond to two qualitatively different wave 
functions? Such dependence on physical conditions in quantum systems may be 
analogous to classical chaotic dynamics. 

Classical chaos is more often seen in discrete systems, which can be described, for 
example, by a discrete map. This was first done by the seminal work of Feigenbaum 
with the aid of iterative maps ([]]] and references therein) where he showed that for 
a large class of iterative maps, exhibiting infinite period-doubling bifurcations, the 
transition to aperiodic orbits can be described by universal constants, governing the 
so-called Feigenbaum route to chaos. Classical chaos corresponds to a disappearance 
of periodic trajectories. As classical chaotic motion is more obvious in time discrete 
systems, it is therefore natural to study the quantum analogy of this phenomenon 
by considering the quantum effects in solids which are naturally discrete in space 
due to their atomic structure. Consequently, we study the atomic lattice taking 
into account a nonlinearity which arises due to electron-phonon interactions. This 
gives rise to the creation of some self-trapped states. Some other systems of non- 
linear lattices arising in film deposition may also be described as coupled nonlinear 
oscillators. 

2 The Hamiltonian And The NSE 

The classical and quantum systems mentioned above may be generally described 
with the use of the Hamiltonian 
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h = Ei^-^r-z|N 4 -££ i^r (i) 

i i i 

where ipi is the wave function of the the self-trapped particle on the i th site, c is 
some parameter and E is the energy eigenvalue. For coupled (nonlinear) oscillators 
the function ip describes a lattice distortion. The wave function, of course, must 
satisfy the conventional normalization condition 

£ I \ 2 = 1 (2) 

i 

which is effectively described when the eigenvalue E is used as a Lagrange multiplier. 
The parameter c may have both positive (for self-trapped quantum states) and 
negative (for nonlinear coupled oscillators) values depending on the system being 
studied. 

The equations which describe the stationary critical points of H are given by 
VH(ip) = 0, where V is the differential operator ie J^. This gives the following 
discretised nonlinear Schrodinger equation (NSE) 

- V'i-i + - ipi+i - c\i)i\ 2 i)i = E^i (3) 

For simplicity, only real solutions are discussed and so we make the substitution 
| ipi | 2 = ipl in the NSE. If we consider the NSE describing nonlinear harmonic oscil- 
lators we have 

- Xi-i + 2xi - x i+ i + cx\ = exi (4) 
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which differs from @ only by the sign of c. With the aid of the transformation 



X n = (-l)> r , 

e = i-E 



(5) 
(6) 



we again get 



- + 2ipi - tpi+i - op* = Eipi (7) 

Thus, we find that all solutions and the classification of such solutions for the 
quantum problem associated with the NSE with positive coupling constant (c > 0) 
are also the solutions (and corresponding classification) of the classical nonlinear 
lattices associated with negative coupling constant (c < 0). 



3 Iterated Maps 

Conversion between Hamiltonian forms and mappings has long been known to be a 
powerful mathematical tool for the theoretical and numerical analysis of dynamical 
systems. Indeed, two-dimensional maps allow representation of a stationary con- 
figuration of a Hamiltonian of the form (|1|) by a trajectory of a dynamical system. 
Therefore, the mathematical situation is identical to that of temporal evolution, 
although the static problem has been formulated in terms of spatial arrangements. 

For a very small number of lattice sites we could attempt to solve the set of 
equations produced for each system. To classify the type of solution obtained we 
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analyse the wave vector structure of this solution. Alternatively, for a larger number 
of sites, we can follow the methods used in the studies of nonlinear dynamical 
systems: we can iterate the NSE ([|), after first representing it in the form of a 2D 
map. That is, we can apply adapted methods from the analysis of chaotic dynamical 
systems and explore the nature of possible solutions to the NSE. To do this we must 
first find suitable orbits. 

The discretised nonlinear Schrodinger equation may be cast to form a first order 
2 — D iterative map in discrete time by introducing the auxiliary variable Zi + \ = 
— ipi (see for comparison 0, [|, |J). This gives 



' Z i+1 \ ( Zi-E^i- C<A 3 ^ 



2 = 0,1,... 



y i>i + Zi+i J 

which can also be considered as a three-term recursion. The total spectrum of 
solutions depending on the initial conditions ipo and Z cannot be obtained in closed 
form. In the early eighties it was common to generate possible solutions of (H) by 
fixing the displacements of two neighbouring particles and then iterating the map, 
however, the physical stability of the system under small perturbations was often not 
considered. In fact, the map technique provides a variety of stationary solutions of 
equation (0). However, the accessibility to all possible solutions is obviously highly 
limited by the parameter range where the map is mathematically stable. Since the 
determinant of the Jacobian of the map equation (||) is strictly one, the map is 
area-preserving and the stability of a cycle ip = {ipi, ip2, ■ ■ ■ , iPn} is guaranteed only 
if the eigenvalues of the Jacobian product along the trajectory are on the unit circle, 
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ie the trace of the Jacobian of the N — th iterate satisfies 



n 2-E- 3c^i 2 -1 



Trn 



< 2 



(9) 



i=l 



1 







V 



/ 



When equation (||) is not satisfied the trajectory {tpi} is unstable with respect to 
the initial conditions. Note that the stability is only local in the sense that for initial 
conditions sufficiently close to the fixed point of the N — th iterate, the correspond- 
ing phase portrait {(^,Zi)} describes the usual elliptic Kolmogorov-Arnold-Moser 
(KAM) curves. Consequently, N— periodic solutions, defined by a fixed point of the 
N — th iterate, depend singularly on the initial conditions and can only be obtained 
by providing the exact initial conditions. 

As an example, consider the fixed point of the evolution of (|8]): inserting (ip*, Z*) = 



The fixed point is only locally stable if —2 < E < 0. 

Note that the stability behaviour does not depend on the parameter c; c can be 
rescaled. The normalisation condition, (Q), can be satisfied by a suitable rescaling 
method with the aid of the free parameter c. If one rescales the wave function %pi by 
/3 and the coupling constant c by 1//5 2 the nonlinear Schrodinger equation remains 
unchanged. In other words, the value ipi and the value c are related by the similarity 
transformation: ^ — > fiipi and c — > c/ (3 2 Therefore we use the scaling parameter 
j3 to satisfy the normalisation condition. The physical value of the (electron-phonon) 
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2 + 2E\ < 2 



(10) 



coupling constant is then C = cJ2i4>i- 

The next step is to fix the values E and C and then iterate the map. Provided 
that for arbitrary values of c a non-divergent solution ipi exists, it will always cor- 
respond to a physical constant C = cJ2i i>j an d the physical rescaled wave function 

is *i = ^/v£W- 

It is possible that for some value of E and a certain lattice size that the iteration 
procedure is divergent. In the classical theory of chaos this is not allowed because 
the phase space must necessarily be bounded. However, in the quantum problem 
this is not important since the system has a finite size. For any finite system we 
simply exclude divergent trajectories from our consideration because this divergence 
means that at such values of the parameters a solution does not exist. 

Numerical experiments investigating different trajectories of the 2-D map have 
been performed for different values of the parameters. We find that there are essen- 
tially only three types of maps produced by iteration. There can be regular maps 
where only a small number of points in the phase space are visited. There are also 
irregular commensurate maps in the form of closed loops. These loops consist of a 
number of points being visited. As the commensurability decreases the points vis- 
ited in the phase space become more and more dense and an elliptic orbit is mapped. 
Finally, we have irregular incommensurate maps where the previous closed loops are 
dispersed (to the stochastic) and a large number of points in the phase space are 
visited. 

These simulations indicate that in this system there are three qualitatively differ- 
ent types of trajectories which depend on the values of the parameters: 1) periodic, 
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2) quasiperiodic and 3) chaotic. It is therefore possible to speculate that the regular 
commensurate solutions may be transformed into the quasiperiodic or chaotic type 
structures by a change of initial conditions or parameters. The type of structure 
which arises in place of the regular one depends on how strong the incommensurabil- 
ity is. For example consider a period 2 solution. If there is an even number number 
of sites then this solution is commensurate with the lattice. However, if there is 
an odd number of sites then there is one site too many for the period 2 solution 
to exist. In this case the incommensurability which is caused by the appearance 
of this site is distributed throughout the lattice but is mostly seen near the extra 
site. If there is a large number of lattice sites then the incommensurability caused 
by the extra site is weak and so the possible solution arising (instead of the regular 
one) will be a quasiperiodic one. The period 2 structure will be slightly deformed. 
This quasiperiodic solution (with respect to period 2) will at the same time also be 
a regular period iV solution if we are considering PBC. This is true for any other 
quasiperiodic or chaotic solution of NSE. 

The incommensurability may be stronger in other cases. For example, consider 
a period m structure on an iV site lattice (m < N) where N is prime. If k such pe- 
riod m structures exist on this lattice then the number of sites can be decomposed 
as N = km + k\ with fci ~ | is a number of the order of k. Then the incom- 
mensurability arising from the extra k\ sites will deform the period m solution. If 
this incommensurability is sufficiently strong then we would expect to have chaotic 
solutions arise instead of the regular period m solutions. 

Below we show that these qualitative arguments are indeed valid. To prove that 
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the proposed classification does exist let us try to find all possible solutions of the 
NSE for a finite system directly applying analytic and numeric methods. 



4 Exact Solutions in the limit c — ► oo 

The above NSE (^) has exact solutions in the limit c —>■ oo for a iV-site lattice, 
which have been published separately || . The associated energy eigenvalue of these 
exact solutions takes the form (for c > 0) 

E= 2m ± U-c 
n 

where the eigenvalue E corresponds to a wave function localized with near equal 
probability on n sites which are separated into m groups (spots). Between the spots 
the wave function is vanishing while inside these spots the wave function changes 
sign a total of I times, although its amplitude is (nearly) the same. Note that in each 
localizing quantum well the wave function is localized with equal probability, ipl = 
1/n for all numbers k associated with localizing quantum wells. This probability 
does not depend on the number of spots, m, the localization area is separated into. 
The number n may take any integer value n = 1, 2, 3, . . . , N. The number m may 
take any integer value satisfying m < n. Similarly, the number I satisfies I < n. 

Eq. (|TTp , obtained in the limit c — * oo, has been compared with the exact and 
numerical solutions for systems consisting of small number of sites. In all these 
cases for nearly all values of c (except small regions of critical values where the self- 
trapped solutions originate) there is perfect agreement with the derived formula (pTTj) . 
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However, in contrast with this perfect agreement between eigenvalues, a decrease in 
the value of c leads to a noticeable deviation in the wave functions (eigenvectors) 
from those obtained in the limit c — > oo. The first order corrections to the presented 
eigenvectors obtained with the use of perturbation theory is of the order of y/n/ c ie 
0(1/ c). When the coupling constant c is not very large the wave functions of some 
states will have interesting incommensurate and chaotic structures. Thus, from 
comparison with numerical results and with perturbation theory we conclude, that 
even though the spectrum of the NSE for a system with a finite, arbitrary number, 
N, of sites is well described by equation (fllD for c 3> 1, the shape of the appropriate 
wave function for smaller values of c may have only qualitative features of the wave 
function obtained in the limit c — ► 00. As c decreases the localization spots smear 



out. Since the spectrum, eq.(|TT|), is associated with a local localization pattern, it 
is universal and does not depend on the boundary conditions. 

Although the structure of the wave function for c ^> 1 does not strongly agree 
with the structure of the wave function in the limit c — > 00 we find that the wave 
function may be approximated reasonably well by an exponential function. This is 
only valid if the peaks in the wave function structure are separated sufficiently so 
that there is little or no interaction between the tails. For lattice sites sufficiently 
far away from the localised site we can assume that in ([?]) the value of the wave 
function is small so that the cubic nonlinear term vanishes. Then, in the continuum 



limit, we obtain a second order linear differential equation which may be expressed 



as 
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Eip(x) 



(12) 



<9x 2 



( HID has the asymptotic solution 



ip = Aexp(—y/—Ex) 



(13) 



where A is some parameter which from the normalisation condition goes as A ~ 
(— E) 1 ^. fljjD may only be used to describe the behaviour of the wave function 
sufficiently far away from the local maxima of the wave function because near the 
maxima the cubic nonlinear term cannot be neglected. Hence, away from the peak 
the behaviour of the wave function may be described as an exponential decay. 

Each of the eigenvalues E correspond to energies created due to a wave function 
localising in quantum wells. Therefore, it is interesting to estimate the Hamiltonian, 
H, which corresponds to the appropriate structures. With the use of the method 
described above the following expression for the Hamiltonian, H, is obtained 



5 Perturbation Theory 

If a small perturbation in the wave function (and energy eigenvalue) is considered, 
with c ^> 1, then a correction to the wave function is revealed. 



H 



c 



(14) 



2n 



Let 
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l/j i= Pi + Xi 



(15) 



and 



E = Eq + E 



(16) 



where pi is the zero approximation wave function, Xi is a small perturbation in the 
wave function, E is the zero approximation energy eigenvalue and E\ is a small 
perturbation in the energy eigenvalue. Substituting this into (^) and ignoring the 
smaller terms gives 

-p i „ 1 -x i „ 1 + 2p i + 2x i -p i+1 -x i+1 -cpf-3cp^x i « E pi + E Q Xi + EiPi + EiXi (17) 
In matrix form the system of N equations becomes 



where F and X are vectors with iV components (N is the number of lattice sites). 
The components of F are Fj = E Q pi + cp\ + — 2pi + Pi+±, the components of 
X are Xi = Xi ie X is the wave vector of first order corrections (~ 1/c) to the zero 
approximation wave function. Finally, (for 3 or more sites) T is a tridiagonal iV x 
iV matrix with diagonal elements = 2 — Eq — 3cp^, elements on either side of the 
diagonal T^-i, T^i + i = —1 and top-right and bottom-left elements T^jy, Tjv.i = — 1 
for periodic boundary conditions (PBC) . The terms involving E\ have been ignored 



TX 



F 



(18) 
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for the meanwhile because \Eq\ ^> \Ei\. Eq, which is given by (|TTD, is in excellent 
agreement with the actual value of the energy eigenvalue for nearly all values of c. 

Premultiplying equation (18) by T _1 gives the first order corrections to the wave 
vector: 



X = T _1 F 



(19) 



Now that we have the first order perturbation corrections to the wave function, 



X{, we can substitute these in (|T7| ) to calculate the first order correction to the energy 
eigenvalue, E. Assuming PBC and summing ( |TTD over all i (neglecting the nonlinear 
terms in X{ as !■ C we obtain the first order correction to the energy eigenvalue 



E 1 



-c Ei p\ - 3c Ei PiXi - E Ei Pi - E Ei Xi 



(20) 



EiPi "t" Ei %i 

To see how the perturbation corrections compare to actual solutions we compare 
these corrections with exact solutions of (0). For a 2-site system with PBC we know 
that the exact solutions correspond to the eigenvalues E = —c/2, E = (8 — c)/2 and 



E = 2 — c. The corresponding eigenvectors are 
for E = -c/2 



i 

V2 



1 



for E = (8 - c)/2 



v/2 
1 

~V2 



14 



for E = 2 - c 




V2 



where a = yl — ^| and lies in the range < a < 1. This solution only exists for 
c > 4. The above three eigenvalues (obtained analytically) are precisely the same 



Note that for E = — c/2 and E = (8 — c)/2 we have no spots ie m = 0. There are 
no sites where the wave function is zero. If Open Boundary Conditions (OBC) were 
considered (instead of PBC) then, due to the nature of the boundary conditions, 

For the solutions which correspond to the eigenenergies E = —c/2 and E = 
(8 — c)/2 we know that there are no corrections to the wave function and energy 
eigenvalue. Thus, we expect that x\ = and E\ = for both of these solutions. 
Substituting for E and Pi in ( PI) , (pOj) we find that this is indeed the case. 

Following the same procedure for the symmetry breaking solution gives no first 
order correction to the localised component of the wave function (but does give a 
correction of — 2/c 2 for c> 1) and a first order correction of 2/c to the vanishing 
component. To compare these results with the exact solutions we must expand the 
surd (in inverse powers of c) for both the localised and vanishing components of the 
wave function. Doing so gives 



as those given by flTT| ) if the appropriate wave function structures are substituted. 




1 



2 



10 



84 



(21) 



V2 
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which agrees completely with the results from our perturbation analysis. Note that 
in the limit c — > oo the corrections are all zero and the zero approximation solutions 
are obtained. Although there is excellent agreement in the wave functions, the same 
is not true regarding the energy eigenvalue. We know that the eigenvalue E = 2 — c 
is exact and consequently there should be no corrections to E . However, we obtain 
the correction E\ ps — 4/c. This is because we did not take into consideration all 
orders of correction. We have only used first order corrections and neglected lower 
order corrections. Clearly, lower order corrections do exist and these too must be 
considered for complete agreement. Thus, the correction in the energy eigenvalue, 
Ei, is zero to an order of c° ~ 1. 

Once we have obtained the first order corrections to the wave function and 
eigenenergy we can then investigate higher order corrections. From the 2-site PBC 
example discussed above we expect that the next correction to the wave function 
will, in general, be of order ~ 1/c 2 . 

6 Iteration Procedure 

We can modify the above procedure by finding the first order correction to the zero 
approximation wave function obtained in the limit c — > oo as detailed above. We 
then find the next order correction to this new wave function. That is we repeat the 
above procedure and determine the next order correction to the 'corrected' wave 
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function. This iterative procedure (which is reminiscent of the Newton-Raphson 
procedure) is repeated until the correction to the wave function becomes negligible. 
This means that the value of the wave function converges to some limit. This limit 
is an exact, numerical solution of the discretised NSE ([?]). 

This exact solution is obtained by following an analogous procedure to that used 
to obtain (103). If ipi(m) is our wave function iterated m times then we again have 

TX = F (23) 

where F and X are again vectors with N components associated with the number 
of lattice sites. The components of F are now F% = Eoipi(m) + cipf(m) + ^_i(m) — 
2ipi(m) +ip i+ i(m) and the components of X are Xi = Xi the (m + l)th iterative cor- 
rections to the wave function iterated to m th order. Finally, T is again a tridiagonal 
N x N matrix with diagonal elements = 2 — Eq — 3c^ 2 (m). 

Premultiplying equation (p3| ) by T _1 gives the (m + l)th iteration corrections 
to the wave vector: 

X = T- X F (24) 
Alternatively, we may write the linear approximation to the wave function as 

Tp(m + 1) = ip(m) + 5ip(m) (25) 

and analogously 
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E(m + 1) = E(m) + 5E{m) (26) 
so that the counterpart of (|T7|) after m iterations may be written in the form 



Vif(*(m + 1)) - V#(*(m)) « V 2 #(*(m))(*(m + 1) - *(m)) (27) 



Then ( f23|) can be expressed in the form 



V 2 #(*(m))5*(m) = V#(*(m)) (28) 
with F = Vi?(*(m)), T = V 2 #(*(m)) and X = 5^{m) so that 

*(m + 1) = *(m) - [v 2 #(*(m))] ~* Vif(*(m)) (29) 

This process is, however, not always convergent. The wave vector may change 
structure when other eigenenergies have a similar value for a particular c. For 
example, for a 3 site system with PBC two of the eigenvalues are — c/3 and (6 — c)/2. 
Both of these have the same value for c = 18 but they have very different structures 
of the wave function. It is possible for the initial wave vector structure corresponding 
to one energy eigenvalue to change during the iteration process to another structure 
which corresponds to a different energy eigenvalue. However, this new structure 
does also exist and is also a solution of the NSE (^). 

Analogous to (^), the energy eigenvalue (for PBC) is given by 



E(m) = (30) 
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and the value of E(m) should converge towards E for most values of c 3> 1 except 
for small regions of c where the localised states originate. Even in these sensitive 
regions the value of E(m) should be in fairly close agreement with Eq. If there is 
a drastic change in the value of E(m) then this indicates that the wave function 
structure has changed and the iteration process will converge to some other limit. 

Let us apply this method to the above NSE for finite c taking the zero approx- 
imation solutions, which are exact in the limit c — > oo, as the starting condition. 
The results obtained by the Iteration Procedure can be compared with the results 
obtained by the use of 2D iterated maps. For such a comparison, however, we have 
to introduce the new phase space function Zi = tp i+ i — ipi for the Iteration Proce- 
dure and plot Zi(ipi) against ipi. This then allows us to compare the wave function 
structures obtained with the use of the Iteration Procedure and the trajectories of 
the 2D iterated maps discussed in Section 3. We find that the results from the two 
different methods are consistent. That is the pictures representing the phase space 
of our quantum system correspond to different periodic or chaotic trajectories of the 
Poincare Map. 

7 Results 

In the previous section we have reformulated our problem of the solution of a discrete 
nonlinear set of equations, (|7|), to an iterative problem. This gives possible criteria 
for the classification of the solutions obtained. Following the Iteration Procedure 
we find that there exist different types of structures of the wave function. We find 
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that at some fixed initial conditions for E > there are the usual plane wave 
Bloch solutions (analogous to the motion of a free electron). The corresponding 
phase portrait has the form of a closed (elliptic) curve. However, nonlinear localized 
solutions arise for E < 0. Our results indicate that there are three qualitatively 
different types of structure of the wave vector for negative values of E: periodic, 
quasiperiodic and (deterministic) chaotic. The same structures can also be obtained 
by numerically solving the NSE for a finite size system. Of course for finite systems 
quasiperiodic and chaotic structures are not well defined because these are, strictly 
speaking, well defined only for systems of infinite size. However, we may still indicate 
analogous features, for example, with the aid of phase portraits adopted for finite size 
systems which are possibly equivalent to commensurate, regular incommensurate 
and irregular incommensurate structures of the wave function on the lattice. 

7.1 Periodic Trajectories 

Consider the initial structure where the normalised wave function is located on only 
4 equidistant sites, ie there are only 4 sites where the wave function is not zero. 
These 4 sites are equally spaced apart on a lattice consisting of a total of 32 sites 
with PBC. Each of the localised sites is separated by 7 sites where the wave function 
is zero. The corresponding energy eigenvalue is given by E = (8 — c)/4. Taking 
the coupling constant to be c = 10 and applying the Iteration Procedure to this 
initial structure of the normalised wave function we obtain the solution presented in 
Figure la. This shows a periodic behaviour of the wave function through the lattice. 
The corresponding phase portrait of the quantum state, Figure lb, indicates that 
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this initial structure gives rise to a regular commensurate period 8 structure. The 
phase portrait consists of only 8 points as the wave function oscillates regularly. 
This solution represents a generalisation of the period two solution for our two site 
system which we obtained analytically. Here the wave function will occupy each 8th 
site with a greater probability than the other 7 sites. Since the lattice deformation 
is proportional to ipf, each 8th site will be deformed by a greater amount. Such a 
regular deformation will give rise to a superlattice of deformations or a superlattice 
of solitons. As c is increased the oscillations in the value of the wave function 
will grow. The width of each peak decreases and the height increases. The wave 
function is a maximum at every 8th site. In the limit c — > oo we will obtain our 
initial structure of the wave function being completely localised on each 8th site and 
zero everywhere else. 

7.2 Quasiperiodic Trajectories 

If the initial structure is altered slightly so that the localised sites are not equidistant 
then we see a deviation from the regular periodic behaviour of the wave function. 
Figure 2a depicts the wave function, obtained via the Iteration Procedure, which 
corresponds to an initial condition for the Iteration Procedure consisting of 4 lo- 
calised sites on a 32 site lattice with PBC. However, this time the localised sites 
are separated by either 5 or 9 lattice sites. The energy eigenvalue is still given by 
E = (8 — c)/4. For c = 12 we see some competition develop between the different lo- 
calised sites. They begin to interfere with each other and give an irregular structure 
which has competing periods but is, nevertheless, commensurate with the lattice. 
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These periods deviate slightly from one another. This is illustrated in Figure 2b 
where we see a dispersion of the points in the phase space. The previous regular 
period 8 structure no longer exists on the lattice. Two different period 8 structures 
arise which correspond to two different sets of points being visited in the phase space. 
The deviation in the period 8 structures indicate the origins of quasiperiodicty or 
even chaoticity within the system due to a strong interaction of different spots or 
different solitons. Specifically, it is easy to imagine a larger lattice with a greater 
number of competing structures. The phase space of the quantum system will then 
become more dense and take the form of a closed loop which would be indicative of 
a quasiperiodic solution. If this closed loop is dispersed then this would indicate a 
stochastic structure which is reminiscent of classical chaotic motion. 

If the value of the energy eigenvalue becomes less negative, ie closer to zero 
but still negative (for example c = 11) then the interaction between the spots 
intensifies. Similarly, if the value of c is increased then the interaction between the 
spots is reduced. In the limit c^oowe will obtain our initial structure of the wave 
function being completely localised on 4 sites and zero everywhere else. Note that in 
this limit the corresponding phase space will consist of just 3 points: one at the origin 
(ipi, Zi) = (0, 0), one at (ipi, Z^ = (0, 1/y/n) and one at (ip u Z t ) = (1/y/n, -1/y/n). 

Figure 3a corresponds to the final state of a wave function which had the initial 
condition for the Iteration Procedure of the wave function being localised on 10 
lattice sites which are all grouped together in one spot on a 100 site lattice with PBC. 
The energy eigenvalue is given by E = (2 — c)/10. For c = 4 the interaction between 
the localised sites gives rise to the wave function being localised in two discrete 
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spots. Note that the value of the wave function is always positive. The phase 
space, Figure 3b, indicates a possible quasiperiodic (or regular incommensurate) 
structure. The phase space of this system lies on a closed loop. The points which 
belong to this phase space appear to be paired. This could be an indication of small 
deviations in the structures of the wave function of the two spots. However there is 
a dispersion of the closed loop trajectory at the origin. This dispersion indicates a 
transition from a quasiperiodic structure to a weakly chaotic regime. The dispersion 
is due to a commensurate wave function structure beginning to transform into an 
incommensurate structure. 

Figure 4a shows the final wave function structure, obtained via the Iteration 
Procedure, which had the initial condition of the wave function being localised on 
12 lattice sites in 12 different spots on a 100 site lattice with PBC. The number of 
empty lattice sites separating the localised sites varies cyclically. There were 6, 7 
and 9 empty lattice sites between different localised sites. This pattern was repeated 
systematically every 25 sites. The energy eigenvalue is given by E = (24 — c)/12. 
For c = 29 we obtain a wave function pattern, Figure 4a, which at first seems to 
have no order. Figure 4b shows the underlying pattern inherent in the wave function 
structure. The phase portrait of this quantum system, Figure 4c, indicates a period 
25 structure. This is due to the structure being repeated every 25 lattice sites. From 
Figure 4b we see that there are actually three peaks of the wave function which occur 
over 25 lattice sites. That is the loop plotted in Figure 4c must actually include 3 
trajectories which span the phase space of the phase portrait. Thus, in Figure 4c 
we see something akin to a regular incommensurate structure of the wave function. 
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7.3 Chaotic Trajectories 

Figure 5a shows the behaviour of the wave function on the lattice sites as c is 
increased to c = 31 for the previous initial condition. The energy eigenvalue 
E = (24 — c)/12 becomes more negative. The wave function structure has split 
into three regions; top, middle and bottom. Figure 5b shows that some type of 
symmetry has been imposed on the wave function. The structure is repeated every 
25 sites as expected. In each set of 25 lattice sites the wave function has three local 
maxima. Two of these maxima consist of just one lattice site but one local maxima 
consists of two lattice sites. This is due to the strong interaction of the spots. The 
corresponding phase space, Figure 5c, is in the form of another dispersed closed 
loop. 

As c is further increased to c = 32 Figure 6a shows the wave function structure 
which has now further separated into very distinct groups. From figure 6b one 
can see that now the interactions between the spots has decreased. The peaks are 
narrower and higher. Figure 6c represents the phase space of this structure which 
is comprised of sets of three points. This chaotic structure is due to the interaction 
of different localisation spots. The appearance of the stochastic layer demonstrates 
a chaotic behaviour and indicates the effects of irregular incommensurability of the 
soliton lattice with the original atomic lattice. Upon further inspection of Figure 6a, 
which shows the value of the wave function through out the lattice, we see that the 
different horizontal regions in Figure 6a actually correspond to the different points 
in the orbit of the phase portrait, Figure 6c. This orbit has 8 sets of three points 
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and on closer inspection of Figure 6a we see that Figure 6a is actually organised 
into 8 horizontal sets; one at the top, two in the middle and five at the bottom. 

Figure 7a shows the behaviour of the wave function on the lattice sites as c is 
increased to c = 36. The energy eigenvalue E = (24 — c)/12 becomes more negative 
and the interaction between the spots decreases. The wave function is again divided 
into three distinct regions at the top, middle and bottom. These distinct regions 
correspond to the dispersed points in Figure 7b. The dispersion of these points has 
now reduced and they are beginning to converge towards different limits. This in- 
dicates that the incommensurate structure is becoming a commensurate structure. 
The wave function pattern in Figure 7a is separated into different sets which cor- 
respond to the different convergence limits in Figure 7b. In the limit c — > oo the 
phase space will converge to the points (ipi, Zj) = (0,0), (ipi, Zi) = (0,1/y/n) and 

(^,Z i ) = (l/v^,-l/ v ^)- 

As c is further increased the wave function structure now begins to resemble the 
initial condition where the wave function was localised on 12 sites in 12 spots. Note 
that this initial structure can only exist in the limit c — > oo. Figure 8a shows the 
wave function of the lattice sites for c = 84. The energy eigenvalue E = (24 — c)/12 
becomes even more negative. The value of the wave function now has three specific 
domains; at the top, near the bottom and at the bottom. These domains correspond 
to the different states of the system as depicted in Figure 8b. The wave function 
structure now indicates regular commensurability. In actual fact, careful study of 
figure 8b indicates a period 8 structure; there are three points closely located at the 
origin. If the resolution of Figure 8b was better it would actually indicate a period 
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8 structure. Note that the sets of three points which caused the dispersion of the 
orbit (see Figure 6c) have converged as c was increased and now there is only one 
discrete point where previously there were three. This dispersion of the points for 
relatively small c could possibly be due to the slightly different corrections of each 
order to the wave function. This is indicated by the 18 possible first order corrections 
(see appendix Al) depending on the structure of the wave function. The first order 
correction for different sites are not all the same. These differences can cause a 
variation in the pattern of the wave function which would result in a dispersion of 
the orbit of a trajectory in the phase space. It is easy to extrapolate this result for 
different orders of correction. For each order of correction to the wave function there 
are an increasing number of possible corrections. More and more different corrections 
appear which decrease for large c and in the limit c — > oo these corrections vanish. 
This then could be the origins of the irregular incommensurability present in our 
quantum system. 



8 Summary 

Thus, we obtain that in the limit c — > oo there arises a degenerate set of solutions 
associated with different localisation patterns of the wave function which consist of 
empty lattice sites (where the wave function is vanishing) and lattice sites where 
the wave function is localised. These localisation patterns, which can be viewed 
as soliton type structures, have many different configurations which correspond to 
the same eigenvalue of the NSE. However, when the value of c is not infinite this 
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degeneracy is broken. Different localisation spots within the pattern start to interfere 
with each other and modify the behaviour of the wave function. The interaction 
between the solitons is determined by the reciprocal of the modulus of an eigenvalue 
of the NSE, 1 / \J\E\, which governs the radius of the soliton peaks. 

The interaction between the solitons depends essentially on the radius of the 
solitons and their separation. If the radius of the soliton peaks is much smaller than 
the separation then we have regular, commensurate behaviour of the wave function 
on the lattice. The interaction between the solitons is very weak. As the radius of 
the solitons increases (1/\J\E\ increases) and the separation remains constant the 
interaction between the solitons increases and the behaviour of the wave function 
on the lattice becomes less regular. Finally, when the radius of the soliton peaks is 
about the same as the separation of the solitons we have strong interactions between 
the solitons. This results in strong incommensurability of the behaviour of the wave 
function on the lattice. 

We find that the interference between the spots can give rise to three qualitatively 
different structures: periodic, quasiperiodic and chaotic. To indicate such structures 
in the quantum state the methods used in the studies of classical chaos were applied 
to the quantum system. A Poincare map for the quantum state was built up. This 
consists of the amplitudes of the wave functions and of residues of these amplitudes 
associated with neighbouring sites. That is it is a projection of the Hilbert 'phase 
space of our quantum system' in the plane, that is, the set {Zi,ipi}. 

In the regular periodic and quasiperiodic solutions the wave function amplitudes 
replicate with some period equal to some integer number of lattice constants or 
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creates aperiodic structures, respectively. However, there is also the appearance of 
structures analogous to the those arising in classical chaos. The destruction of the 
periodic and quasiperiodic orbits, which has been ascribed to the creation of chaotic 
structure in the wave function, is also exhibited by the system. That is the creation 
of a structure which has no definite period has been found. 
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Figure Captions 



Figure 1 The final (converged) behaviour of the wave function on the lattice (a) 
and the corresponding phase portrait (b) obtained by Iteration Procedure for a 32 
site lattice with Periodic Boundary Conditions (PBC). The initial conditions were 
c = 10 and E = —0.5 with the initial, normalised wave function value of ip = 1/2 at 
only 4 equidistant sites. Between these 'localised' sites the value of the wave function 
was zero. Figure la shows a regular (periodic) behaviour of the wave function. The 
pattern in Figure lb indicates that only 8 points are visited all of which lie on a 
closed loop. Therefore, this structure of the wave function corresponds to a regular 
class of behaviour. 

Figure 2 The final behaviour of the wave function on the lattice (a) and the 
corresponding phase portrait (b) obtained via the Iteration Procedure for a 32 site 
lattice with PBC. The initial conditions were c = 12 and E = — 1 with the initial, 
normalised wave function value of ip = 1/2 at only 4 sites. These 4 sites were not 
distributed uniformly throughout the lattice but had either 5 or 9 sites where ip = 
giving an 'irregular' structure. Figure 2a shows the irregular behaviour of the wave 
function. Figure 2b shows the 16 points which are visited in the phase space. The 
points no longer lie on a closed loop and so this arrangement of points is associated 
with an irregular class of behaviour of the wave function. 

Figure 3 The final behaviour of the wave function on the lattice (a) and the 
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corresponding phase portrait (b) obtained via the Iteration Procedure for a 100 
site lattice with PBC. The initial conditions were c = 4 and E = —0.2 with the 
initial, normalised wave function value of ip — 1/y/lO at 10 sites. These 10 sites 
were all gathered together in one spot (group of localising sites). Figure 3a shows 
the behaviour of the wave function after it has converged. Figure 3b shows the 
trajectory of this wave function. Again a closed loop is traced out by the points 
visited in the phase space, however, this time there is some dispersion of points 
at the origin. The loop appears to consist of pairs of points; each point in a pair 
corresponding to the change in the wave function at equivalent locations in either of 
the two separate peaks observed in Figure 3a. This arrangement of points in Figure 
3b indicates some type of irregular structure of the wave function which is somehow 
commensurate with the lattice. 

Figure 4 The final behaviour of the wave function on the lattice (a) and the 
corresponding phase portrait (b) obtained via the Iteration Procedure for a 100 
site lattice with PBC. The initial conditions were c = 29 and E = —0.416- with 
the initial, normalised wave function value of ip = 1/ V^12 at 12 sites. These 12 sites 
were separated into 12 spots which were further separated by either 6, 7 or 9 'empty' 
sites. This structure was repeated every 25 sites. Figure 4a shows the apparantly 
random behaviour of the wave function after it has converged. Figure 4b shows 
the underlying pattern manifest in this structure. The pattern is repeated every 
25 lattice sites as expected. Figure 4c shows the trajectory of this wave function. 
Once again a closed loop is traced out by the 25 points visited in the phase space. 
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However, this loop corresponds to three distinct orbits; one for each of the three 
peaks in the wave function every 25 lattice sites. 

Figure 5 As Figure 4 but with c = 31 and so E = —0.583'. Figure 5a shows the 
apparantly random behaviour of the wave function after it has converged. Figure 
5b shows the underlying pattern manifest in this structure which repeats itself every 
25 lattice sites. Figure 5c shows that the 25 points visited in the phase space no 
longer lie on a closed loop. However, most of these points seem to belong to a set 
of three points. This could again (see Figure 3) correspond to the change in the 
wave function at equivalent locations in the separate peaks observed in Figure 5b. 
This arrangement of points possibly indicates some type of irregular structure of the 
wave function which is incommensurate with the lattice. 

Figure 6 As Figures 4 & 5 but with c = 32 and so E = —0.6'. Figures 6a and 6b 
show the behaviour of the wave function after it has converged and the underlying 
pattern manifest in this structure, respectively. Figure 6c shows the 25 points visited 
in the phase space. It can now be clearly seen that most of these points occur in 
sets of three. As the peaks in the wave function get narrower and higher (Figure 
6b) the dispersion of these points decreases. The behaviour of the wave function 
in Figure 6a seems to be divided into zones at the top, middle and bottom. These 
zones correspond to the sets of different points in Figure 6c. 

Figure 7 As Figures 4, 5 & 6 but with c = 36 and so E — — 1. Figure 7a shows 
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the behaviour of the wave function after it has converged. Figure 7b shows the 25 
points visited in the phase space. The zones in Figure 7a are much more distinct 
and further apart as c increases. 

Figure 8 As Figures 4, 5, 6 & 7 but with c = 84 and so E = —5. Figure 8a 
shows the behaviour of the wave function after it has converged. Figure 8b shows 
the points visited in the phase space. 

Appendix Al 

As the first order corrections to the wave function depends mainly on the value 
of the zero approximation wave function at only three neighbouring sites it is simple 
to obtain algebraic expressions for all first order corrections to the wave function. 
Note that, although the first order correction to the wave function at a particular 
site does have some dependence on the first order corrections to the wave function 
of its neighbouring sites, because <C and x i+ i <C p i+ \ this dependence is 
weak and can, often, be ignored. In QT7D, Pi^i,pi,p i+ i may be either 0, or 
—1/y/n (assuming x^i,x i+ i ~ 0). 

For pi = 0: 

Xi = if both pi_i and p i+ % are 0. 

Xi = r^z^zm ^ either pi_\ is 1/y/n and p i+ % is or vice versa. 
Xi = ( c _2^ 4 ;) if either is —1/y/n and pi + \ is or vice versa. 
Xi = (c-l^-41) i f botn and p i+x are 1/y/n. 
Xi = if either is 1/-Jn and p i+ i is —\jyfn or vice versa. 
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For pi = 1/y/n: 

X i = Ji[cTm+2l) if b ° th aIld are - 

Xi = 2 ^(l+m+2i) ^ either pi_i is 1/y/n and is or vice versa. 
x i = 2%{c+m+2i) if eitn er pi_ x is -1/y/n and p m is or vice versa. 

X i = if b ° th P*- 1 alld ale Vv 7 ™- 

x * = ^{c+m+2i) if either is and p i+ i is -1/y/n or vice versa. 

x * = ^"^"+20 if both and are - Vv 7 ^- 
For = 

x i = - ^(c+m+20 if both p^ and p i+1 are 0. 

x * = - 2^c+m+2i) if either p;_i is 1/y/n and is or vice versa. 
x i = ~ 2vS"a) if either p^ is -1/v 7 ™ and p i+1 is or vice versa. 

X i = ~ vl(c~+m~+20 if b ° th P^ 1 aIld ^ 

Xi = ~ Ji(c+m+2i) if either p;_i is 1/v 7 ™ and p i+1 is -1/y/n or vice versa. 

^ = ^R™+m+2l) if b ° th aIld ale -f/v 7 ^- 
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